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ABSTRACT 


The present work is an attempt to carry out the nonlinear 
analysis of plates and linear analysis of shells. In both the 
cases T transverse shear deformation have been included. The 
nonlinearity is introduced by von-Korman strains. The response is 
obtained using first order shear deformation theory. Solutions 
have been obtained using finite element technique, ft nine noded 
doubly curved isoparametric lagrangian finite element is used to 
d scretize the structure. The results are presented for linear 
bending, nonlinear bending, linear buckling and postbuckling of 
plates. However for shells only linear bending and buckling 
analysis have been conducted. 

It is found that the first order theory (FSDTJ may involve 
about S.5% difference with the higher order theory (HSDTJ for a/h 
= 2 . After a/h = 10 there is hardly any significant difference 
between FSDT and HSDT which requires" lot more computational 
effort- It is also found out that in case of buckling, FSDT theory 


gives resiilts which are close to 

Reddy ’s 

HSDT 

results 

for 

all 

thickness ratios- This suggests 

that one 

need 

not go 

for 

the 


higher order theory in case of the buckling. 

To solve the nonlinear equations in the nonlinear analysi 
an incremental Newton—Raphson method has been developed which i 
more efficient than Picard type of iteration- The theoretica 
results are compared with the experimental results and are foun 
to be in close agreement. The effect of different simply support 



boundary conditions has also been investigated. 

Separate programs written in FORTRAN are developed -for 
symmetric and unsymmetric laminates to reduce the computational 
effort . 



CHAPTER 1 


INTRODUCTION AND LITERATURE REVIEW 

1.1 lnt,roiiuct,Qry Remarks : 

In recent years fibre reinforced composites^ have found 
increasing use in many engineering structures and consumer 
applications, ranging from fuselage of an airplane to the frame of 
a tennis r^uet. This is mainly due to two desirable features of 
reinforced composites ie., a high strength to weight ratio and 
fieKibility in tailoring the anisotropic material properties by 
changing the fibre orientation, stacking sequence and thickness of 
laminae. The increasing use of composites as a structural element 
has generated considerable interest in the analysis of composite 
plates and shells. 

Composites can be considered as a material consisting of 
two or more chemically distinct constituents with a distinct 
interface between them. The properties of composites are 
noticeably different from constituent materials. Fibrous 
composites consist of essentially high strength fibres embedded in 
a matrix. A lamina is a thin sheet with fibers, oriented in 
certain directions, forming the main load bearing members. The 
matrix keeps fibres in position and transfers the stress to fibres 
through bonding. The lamina can be considered to behave as a two 
dimensional orthotropic material. 

The main advantage of composites lies in tailoring the 
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material properties to suit the applied loads. The disadvantage, 
from analytical point of view is the existence of coupling between 
bending and stretching. Due to this the governing partial 
differential equations for laminae become coupled in displacements 
and rotations. This boundary value problem is difficult to solve. 
Very few exact solutions can be obtained for limited boundary 
conditions. This is a marked departure from conventional isotropic 
materials. 

Due to the high ratio of elastic modulus to shear modulus 
(of the order of 25 to 40, instead of 2-6 for a typical isotropic 
material), shear deformation can be significant in composite 
materials as compared to isotropic materials. The classical plate 
theory which does not take into account the transverse shear 
deformation, under predicts the deflections and over predicts the 
buckling loads and the natural frequencies, when the ratio of 
lateral dimension to the thickness is less than 50. 

For small deflections ie., the deflections which are small 
compared to thickness, the linear theory gives sufficiently 
accurate results. Plates, subjected to lateral loads and/or in 
plane loads, fail at large deflections. For large deflections the 
lateral deflection is accompanied with mid plane stretching. This 
increases the load carrying capacity of the plate. Associated 
increase in load depends on the in plane boundary conditions 
imposed on the edges. This behavior is extremely important in 
structural weight minimization in aerospace structures- A plate 
starts to buckle when the critical load is reached. The maximum 
stress at the instant of buckling is much lesser than the yield 
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stress of the material. Hence the plate can take more load. As the 
in plane load is increased* the deflections increase at a faster 
rate. Now the deflections are no more comparable with the 
thickness, violating small deflection theory. Therefore large 
deflection theory should be considered to determine accurate load 
carrying capacity of laminated plates. 

The rapid development of computing technology has 
completely changed the out look of structural analyst. Of all 
approximate methods to analyze the structures, involving complex 
geometries and variety of boundary conditions, finite element 
method (FEM) is found to be most powerful. FEM, in structural 
analysis, is a technique that discretizes a structure in to a set 
of simple sub domains called finite elements. Over each element 
the solution of the governing differential equation is 
approximated by a linear combination of undetermined parameters 
and pre selected algebraic polynomials. The undetermined 
parameters are determined such that the variational statement 
obtained from the governing differential equations is satisfied in 
each element. The satisfaction of the variational statement 
results in a set of algebraic relations among the undetermined 
parameters of an element. Since the domain is divided into a set 
of elements, they are put together into their original position, 
called assembly of elements, using continuity at element 
interface. The finite element method typically entails the 
solution of a very large number of equations for the nodal values 
of the function being sought. Thus the method requires so much 
computation that it is practical only if calculations are carried 
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on a computer. 

In the present work FEM is used to carry out the linear and 
non linear analysis of composite laminates including transverse 
shear deformation. 

1.2 Literature Review : 


The pertinent literature survey is presented here under the 
following three headings. 

(i) Studies using small deflection theory of plates. 

(ii) Studies using large deflection theory of plates. 

(iii) Studies devoted to small deflection theory of shells. 

1.2 Ci]> Studies Using Small Deflection Theory of Plates i 


A number of theories have been proposed for shear 
deformable laminated plates. Stavsky cm extended the classical 
plate theory for isotropic materials to laminated plates. This 
theory has been generalized to laminated anisotropic plates by 
Yang, Norris and Stavsky LZl and is called as YNS - theory- This 
is generalization of the Mindlin plate theory which includes 
transverse shear deformation and rotary inertia effects- In this 
theory the normal to the mid plane before deformation remains 
straight after deformation but not necessarily normal to the mid 
plane. To make shear strains zero on top and bottom faces of the 
plate, a shear correction factor is used which depends upon the 
ply lay~up and fibre orientation, etc- It has been shown that 



C3-53 YNS theory predicts the transverse deflections and natural 
frequencies for laminated composite plates with reasonable 
accuracy. Whitney and Pagano C63 presented some close form 
solutions for cross ply and angle ply rectangular plates under 
sinusoidal loading. Sinha and Rath £73 considered both vibration 
and buckling for the same type of plates. They investigated the 
influence of shear deformation, panel curvature and aspect ratio 
on bending stretching coupling. They used Navier type solution to 
reduce the governing partial differential equations to algebraic 
equations. 

The finite element analysis of laminated composite plates 
began with Pryor and Barker C83 who employed twenty eight degree 
of freedom (dof) shear deformable rectangular element. Pica C93 has 
compared the results of linear, serendipity, lagrangian and 
heterosis elements with analytical thin plate solutions for 
square, circular and elliptical isotropic plates based on YNS 
theory. He suggested the use of heterosis elements for rectangular 
and curved element meshes. Reddy C103 developed C° penalty plate 
bending element based on YNS theory. The element consists of five 
dof at each node- three displacements and two rotations at the mid 
plane of the plate. Reddy C113 presented comparison of the finite 
element solutions using penalty method and exact solutions for 
anisotropic plates. He used reduced integration to correct the 
shear stiffness. He concluded that higher order elements are not 
affected by locking. He showed that nine noded lagrangian element 
is not much affected by the integration scheme but eight noded 
element is profoundly affected- The consideration of syimetries in 
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unsymme'tri cal laminaled composite plates was investigated by Reddy 
11123. He presented the boundary conditions to be used in -full and 
quarter plate models of cross ply and angle ply laminated plates. 
Hass and Lee C133 used a nine noded degenerate shell element based 
on modified Hellinger Reissiner principle. Functions were assumed 
independently for in plane and transverse shear strains. 

The finite element formulation for the isotropic plate 
buckling, using classical theory, is available in standard texts 
C14,153. Gajbir Singh and Sadasiva Rao C223 used an eight noded 
Mindlin element for the stability analysis of thick angle ply 
laminates. The effect of fibre orientation, material properties, 
layering and boundary conditions on buckling loads were studied. 

Further development in laminate analysis have been higher 
order theories. Phan and Reddy C163 presented a higher order 
theory where in plane displacements are expanded as cubic function 
of the thickness coordinate while the transverse deflection is 
kept constant over the thickness of the plate. They presented 
closed form solutions for the bending of symmetric cross-ply 
laminates. Recently Phan and Reddy C173 using the higher order 
theory presented closed form solutions for the stability of 
isotropic, orthotropic and laminated plates. 

1.2 Cil3 S^udlms Itelng Largm Def'ormat.ion Thoory of Pla-boB i 


Much of research in the laminated composite plates is 
limited to linear problems. This is due to the complexity of non 
linear partial differential equations associated with large 
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deflection theory of contposite plates. Approximate solutions for 
composite plates using large deflection theory based on Von Herman 
theory of isotropic plates have been attempted by C23-S43. Chandra 
and Raju tI233 and Chia and Prabhakara C243 employed the gelerkian 
method to reduce the governing differential equations to an 
ordinary differential equation, in time for the mode shape. The 
perturbation technique was used to solve the resulting equations. 
They analyzed bending and free vibration of composite laminates. 

Using FEM the total stiffness matrix can be written as the 
sum of linear stiffness matrix and non linear stiffness matrix. 
The non linear stiffness matrix is called as incremental stiffness 
matrix. Two notations have emerged for the incremental matrices 
formulated in total lagrangian system. In total lagrangian system 
all the quantities of interest are expressed as functions of 
rectangular co-ordinates of the body in its initial or undeformed 
state. The first notation is due to Mallet and Marcal C253 and 
will be referred to as ENT notation. In this type of formulation, 
the expressions for total energy, and equilibrium are expressed as 

T T 

n = {^} - 4 [ ] )W - B {«} 

[ [ 3 + 4 [*^3] B ' B 

(-0 

Where K is the linear stiffness matrix and elements of N1 
*- 

and N2 depend linearly and quadratically on the transverse 
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deflection u^. 

The incremental equilibrium equations can be written as 

[ [ « ] + [ N2 ] ] = {ar} (I.2) 

The second notation is due to Zienkiwicz and Nayak C263 

referred to as CBD notation In this type of formulation the 

equilibrium equation written as 

[ [ X]. - [ * [ X].] {r} = {«} 0-3) 

where K is the linear stiffness K is unsymmetric stiffness 
matrix and K is the initial stress matrix, due to large 

€f 

deformation theory. 

The incremental equilibrium equations can be written as 

[ K {q} = {ar} Cl 4) 

where is the tangent stiffness matrix. The various matrices K^, 

K , K and K are derived for first order shear deformation 

MU T 

theory in chapter 3. 

The correlation between the two notations was presented by 
Wood C273. 


Reddy and Chao C283 used mixed finite element for non 
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linear bending of composite plates in total Lagrangian system. 
Putcha and Reddy CS93 used mined finite element for the non linear 
analysis of laminated plates based on higher order deformation 
theory. They used Picard type of iteration to solve the non linear 
algebraic equations with small load increments. 

1.2 Studies Etevoted to Small Deflection Theory of iS^ells : 


The analysis of laminated composite shells began with 
Ambartsumyan C303 . He assumed that individual layers are oriented 
so that principal axes of material symmetry are coincident with 
the principal coordinates of shell reference. Dong, Pister and 
Taylor C313 formulated a theory for laminated shells which is an 
extension of the developed by Stavsky C13. His work can be called 
as Donnel's shallow shell theory. Cheng and Ho 11323 analyzed 
laminated anisotropic shells using Flugge*s shell theory 11333. The 
theory which considers shear deformation and thermal expansion 
through the thickness was considered by Zukas and Vinson C343. 

The finite element analysis of composite shells was 
attempted by Pandya and Natarajan C353 and Shivakumar and 
Krishnamurthy i:363 and Rao K.P. C373 based on Kirchoff's theory. 
Pandya and Natarajan i:3S3 used a three noded one dimensional 
curved finite element having four degrees freedom at each node for 
the analysis of axisymmetric fibre reinforced laminated shells, 
subjected to axisy^etric loading. The results for a cylindrical 
tube subjected to internal pressure and axial tension, are 
compared with the analytical solution. Shivakumar 


and 
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Krishnamur'thy iI3£»Zi used eight do-f eind sixteen do-f ring element •For 
free vibration analysis of laminated cylindrical and conical 
shells. They showed that for the same number of unknowns results 
with sixteen dof ring element are in good agreement with the 
experimental results than eight dof element. 

Reddy [I38II developed shear deformable finite element for 
composite shells based on Sander's theory C19I]. He compared the 
approximate shell theories of Love , Sander, Donnel, Merely and 
Loo. He concluded that the effect of shear deformation on the non 
dimensional deflections is more pronounced with the decrease in 
layers. He showed that there is a small difference in results 
obtained with different shell theories for shallow shells 
<R/h > 100). But for R/h = 10, the difference in the results 
obtained with different shell theories is noticeable. Using Von 
Korman strains and Sander's theory, and incorporating transverse 
shear deformation, Reddy et.al.C401 presented geometric non-linear 
analysis of doubly curved shells. 

1.3 Aim and Skiopa of PrasanL Work i 


The literature survey presented reveals the necessity of 
nonlinearity and shear deformation in the composite plates. There 
are some higher order theories proposed till now to improve the 
solution accuracy for the linear problems. These are 
computationally expensive and involve large number of nodal degree 
of freedom. The analytical methods can be applied to linear 
problems with limited boundary conditions. But they are not 
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tractable for nonlinear problems. 

The scope of the present study can be summarised below 

In Ref. Cim it is stated that higher order elements are 
not much effected by locking and nine noded lagrangian element is 
superior than eight noded element. Therefore the present study is 
an attempt to study the analysis of plates and shells using nine 
noded shear deformable element. The displacement functions require 
C° continuity which leads to easy formulation for thick and thin 
plates. 

To compare the numerical reliability and numerical accuracy 
of the first order shear deformation theory CFSDTII compared to 
higher order theory CHSDT3 for linear bending and buckling 
problems. The FSDT is computationally inexpensive and easy to 
implement compared to HSDT . 

To check the validity of the Von-Korman nonlinear theory, 
the theoretical results are compared with the experimental 
results. 

A pi card type of iteration is used in most of the works in 
solving the nonlinear equations [128,2911. This type of iteration is 
computationally expensive. Therefore incremental Newton-Raphson 
method has been developed. 

In the nonlinear analysis the effect of different in plane 
boundary conditions is significant- The effect is more pronounced 
with bending— stretching coupling. Hence the nonlinear analysis of 
two layer cross-ply and two layer angle-ply with S1,S2 and S3 
boundary conditions given in chapter 2 is studied. 
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1.4 Layout of tho Thesis : 


The thesis has been divided in to five chapters. The first 
chapter is an introduction to the problem with the pertinent 
literature. 

The second chapter deals with the theoretical formulation 
and derivation of the expressions which can be readily used in the 
finite element formulation. THe finite element formulation and the 
solution procedure is presented in chapter 3. Chapter 4 deals with 
the interpretation of the results obtained and is followed by 
conclusions and scope of further work in chapter 5. 

The thesis has Appendix— 1 and Appendix— 11. The figures and 
tables are given in the chapters that follow as and when the need 


arises. 



CHAPTER a 


THEORETICAL FORNftILATION 


a. 1 Introduction s 

This chapter presents a detailed formulation of shear 
deformation theory of shells based on Sanders theory, including 
geometric nonlinearity. The geometric nonlinearity is introduced 
by von-Korman strains. The nonlinear partial differential 
equations for plates and shells cannot be solved analytically for 
complex geometries and various boundary conditions. FEM is the 
most effective numerical method to solve such problems. 

FEM is a piece wise application of variational method. In 
variational method the governing differential is recast in an 
equivalent integral form and a variational formulation is derived. 
This results weaker continuity requirements on the solution in the 
variational problem than in the original equation. 

In solid mechanics problems the variational form of the 
governing equations represents the total potential energy of the 
system. The solution of the differential equation is found by 
minimizing the variational form with respect to unknown variables. 
For monotonic convergence the primary variable should be 
continuous up to one order less than highest order, occurring in 
the variational functional. 

The classical plate requires cubic or higher order 
interpolation functions for approximating transverse deflection. 
But this first order shear deformation theory allows us to use 
linear or higher order interpolation functions for all 
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- This lype of f ormula'tioin is easy to iiaplei&ent and 
applicable for both thin and thick laminates. However, shear 
deformable elements suffer from locking phenomenon. It is 
necessary that the shear strain fields, as the plate thickness 
reduces, should satisfy the Kirchoff’s condition ie.. 



— ^ 

4- w 

> o 

xz 

K 

yK 


and y 

=: ^ 

Hh ^ 

> 0 

yz 

y 



Conventional methods 

using 

enact 

integration 


consistently, introducing spurious constraints that lead to 
locking and severe stress oscillations. Therefore reduced 
integration technique is used to circumvent the locking. 


2. 2 Assun^tlons made In the Theory t 

Theoretical development in the present work is based upon 
the following assumptions. 

(i) The laminate is constructed of an arbitrary number of 
layers of orthotropic sheets bonded together. However the 
orthotropic axis of material symmetry of an individual 
layer need not coincide with the laminated coordinate 
system. 

(ii) The strains are small. 

(iiilThe transverse normal stress and strain are negligible. 

(iv) The line normal to the mid plane after deformation remains 
straight but not necessarily normal to the mid plane. 

(v) There are no body forces. 

Cvi) Stress-strain relations are linear. 



2. 3 Geomet^ry : 


-i-D 


Fig. 2.1 shows a differential element of a shell. Here 
C ^ » ?g» C ) denote curvilinear coordinates such that and 

^gcurves are lines of curvature on the mid surface C = O- The C 
curves are straight lines perpendicular to surface C = O. The 
lines of principal radii of curvature coincide with the coordinate 
lines and the values are denoted fay and Rg. and Xg denote 

cartesian reference frame. 

The linear displacement field is given fay, 

u ( K ) = u^ ( K O ) + C 

V ( f eg, c > = Ug < e^, eg, o > + e ^2 <2-1) 

w ( e^» eg' c ^ ~ “3 * ^1' ^2' ^ ^ 


Where u, v, w are the displacements along <e^T eg' coordinates, 
Ui , Ug, Ug are the displacements of the corresponding point on the 
middle surface and and ^g are the rotations at ^ = O about eg 
andC^^Kis respectively. 

The strain displacement relations including geometric 
nonlinearity can be written as 


■1 

‘2 

'4 

'5 

‘6 


o 
= 1 
ja 

o 

=4 

o 

“5 

& 


K 

e 


+ e K 


( 2 . 2 ) 
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yhere «. and k. are the mid plane strains and curvatures 

Vi. ^ 

respectively. They can be written as C403 
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.C2.3a) 


Where 




1 

R. 


) 


C2-3b) 


The term C is due to Sander's theory which accounts for the 
o • 

condition of zero strain for rigid body motion. 
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2-4 SJtr-ess -Strain Ralatlons tor- a Lamina * 

First order shear deformation theory used in the present 
work is incapable of satisfying the condition of zero transverse 
strain on top and bottom surfaces. Correction factors K^and 
have been used to account for this in stress strain relations. 

Based on the assumption that the transverse stress is 

negligible, the stress strain relations for an orthotropic lamina 
k oriented arbitrarily to the laminate axis are given by 



Eqn.(2.4} can be written in more compact notation as 

Expressions for the elements of CQ3 matrix in terms of engineering 
constants for lamina for a given orientation « are given in 


Appendix— I- 
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Z.S Stress— RtifSul'ban.trS at< t,he Mid Plane of* Ihe Laainat.e s 


The stress resultants , Ng,N^, , Mg, M^, and Qg of 
the laminate with nl number of layers is given by 


( N^, Ng, ) 


Hn 

E 


K-M 




, <>-g, > d C 


K-*i 


1' E 



Hla 

r 

M, ) = 

o 

E 

I ^ ' "^2' "^6 ^ ^ 



r 

) 

E 

I < ^4- ^5 ^ d C 




resultants 

are 

as shown in Fig.E.2. 


(E.5) 


Substituting relations (E.E) in Eqn.2.4 the following 
expressions for stresses are obtained 




1 

®12 

®16 ■ 


-V 


■ °11 

®1E 

«16 ■ 


f k 1 


k. ss 


22 

®26 

- 

a 

'“e 



Qgg 

«26 

< - 


^ "^6 - 


Sym 


®66 - 


a 

L J 


Sym 


®6& - 


k* 

L J 



C2,6) 


& 


4 


-^5 


y = 


Q 

1 44 


®45 

•‘I SsJ 




Rewriting the stress resultants CS.5) in terms of mid plane 
strains and curvatures after integrating through each layer 


^ ■ 



r 1 


r -| 

■ 


1 

*^2 



QA 


1 OB 



& 



L J 

. 

L J 



^2 

•^6 

- 



m M m m 

• 

« ■ ■ 

m « 

m m 

H 

o 





« 

m 

■ •] 




»2 



Sym 

M 

m 

OD 



4 

^ «6^ 




m 

m 




• 4 ' 


Cl ■ 


QS 


Cl) 


(2.7) 


Where GM, QB, QD and OS are extensional, bending - 
extensional coupling, bending, and shear stiffness matrices 
respectively. They are given by. 


r [ QA], [ OB], [ QD] 1 = E J [ Q] ( 1. C, dC 

i = i = 1,3,6 (2.8) 


1^3 


JMl* ^ 

E J [^] <iC 

X.*4 t 

^ 1C 

i = j = 4,5 


(2.9) 
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£.6 Boundary coiuiltlons s 


Since displacement based formulation is used in the present 
study only displacement boundary condition has to be specified. 
One can develop Navier solutions for different types of laminates 
with three types of simply supported boundary conditions. They are 
referred in this thesis as SI, S2 and S3 respectively. 

The functions which satisfy S1 boundary conditions are 



F ^ 

li 

£ 

mas 4 

U €li>& 

mri 

1 

SIH 

/3Xg 

UgCx^ 

F Kg ) 

= £ 
nssi 

E 

m«4 

V SXM 

mr» 


OOS 

/5Xg 

UgCX^ 

FXg) 

= £ 
ri*A 

£ 

m*4 

y £SXH 

mrt 


SXN 

/^Xg 


FKg) 

= £ 
r»=i. 

£ 

X COS 

mr* 

atx , 
1 

SXM 

/5Kg 

^2<«1 

F Kg) 

= E 

ri3ss4 

z 

m*4. 

y SIN 

mri 

1 

COS 

/5Xg 

N^(o, 

K2) = N, 

(aFXg) 

= «2 

(X^fO) 

” '^<“1 

,b) 

M^(o, 

Kg) = 

(aFXg) 

= ”2 

<X^fO) 

= 

Fb) 


iz.soy 


Where a 
shell . 


= — and = •; — f a and b are plane form dimensions of the 

a b 

The expressions which satisfy S2 boundary conditions are 


“I'N 

FXg) 

= z 

E 

wi™* 

U SXH 

mrt 

OIX^ 

1 

cos 

/5x 

“2**1 

FXg) 

* Z 

',n*4 

Z 

»«■* 

V oos 

mri 

1 

mxM 

/Sx 


N^(o,Xg) = N^CarXg) — N^(x^,o> - N^(x^,b) = O 

M^CoFXg) = M^taFXg) = MgCx^Fo) = HgCx^Fb) = 0 
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The remaining expressions are same as given in Eqn.S.8 

The boundary conditions S1, SS and S3 for full shell 
quarter shell for m * n = 1 are as shown in Fig. 2. 3 

For clamped case all the displacements and rotations 
restrained at the edge ie., 

Ujj — Ug = Ug ~ ~ 


and 


are 


0 


( 2 . 11 ) 








CHAPTER 3 


FINITE ELEMENT FORMULATION AND SOLUTION PROCEDURE 

3. 1 Introduction : 


In this chapter a displacement finite element formulation, 
based on total lagrangian for nonlinear analysis is presented. The 
continuum is divided in to number of elements and over each 
element a set of functions are chosen to define the state of 
displacements. The total potential energy is written in terms of 
displacements and using the principle of stationary of total 
potential for an element, we obtain equations which are nonlinear. 
These equations are solved by Newton-Raphson method. 

3. 3 Element s 


The element used in the present work is a shear deformable, 
doubly curved shell element whose projection is an isoparametric 
rectangular element on plane. Over the typical shell 

element A , the displacements u^, Ug, Ug, and Hjt-g are 

interpolated by the expressions of the form 

HN 

I N. < < 

J=-l 


L 


1,3 


(3.1a> 
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MM 

e 





* 

= 2 < ‘fi' 'a' < 


i == 

1,2 

{3.1b> 


J = 1 





Where mn is the 

number of nodes (9), N 

is 

the 

shape 

functions 


J 

e 

e 




associated with 

node j and u , <* ** 

i i 

are the 

unknown nodal 

variables. The 

shape functions are given 

in 

C143 




af. 3 Total Lagrangian Virtual Work Equation s 


The basis of the present formulation is the virtual work 
equation for a continuum written in total lagrangian coordinate 
system under the assumption of small strains. 


I d ^ ■ 13 '. dv 


1 p du.’^ q. dv 

* r 

du.'*' R. dA 

J " ^ 


J i 

J 

a ^ 


* i-1,5 

C3.2) 


Where v is the undeformed volume, p is the mass density, q is the 

body force per unit mskss end R. is the surface tractions per unit 
areay acting over the undeforined area* Here sy»bcil d is used to 
denote first variation. 

The internal virtual work d fct® over a element is given fay, 

T 

d dv C3«3a) 
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hy-E 


d W 


J J 


d C dA 


C3.3b) 




Substituting the strain relations Eqn.(S.S> in to 
Eqn.{3.3b}, we obtain. 




dW 


“II [ <“ '‘V 


+ d( «g + C kg > + 


- h >'2 


d( + <fk'^>o‘.+ d <D'- + d «* O'- d C dx_ 

6 66 44 55 '1S 


(3.4) 


Where h is the thickness of the laminate. After integrating 
through the thickness of the laminate, Eqn-3.4 can be written as 


dW. 


J [” 


^ d «" + dk® + Ng d «2 + Mg dkg + 


*^6 ”6 d + Qg d 4 ]dx^ dx. 


(3-5) 


Defining the mid plane strains, bending strains and shear strains 
in a matrix form as. 
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{‘F 

{‘F 


^2' ^6 


{ ‘i- 

{ *'"■ “V '^6 } 

{ %■ ‘5 } 


(3.6) 


The mid plane strain can be written as sum of linear 
strain and nonlinear strain. The linear strain is due to membrane 
strain and the nonlinear strain is due to bending deformation u„ 


only. 


Where 


{^r 




<3. 7a) 


{-):■ 




^ 1 

* 2 ^ 
3 l 




— 

dx^ R . 

1 1 


2 

J 

49u_ u„ 

r 

1 

^“3 1 ^ 



— 

— 


V Jni. 

2 

L dXg J 

dll dUg 


f ^ I C ^31 

. *“2 . 


1 #x^ J 1 dXgJ 


(3.7b) 


3.3.1 EKpr«!Ssi€3ffii For i 

"i- 


'K 


- b : 


The membrane strain can be written in a matrix form in 

^1* 
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terms of unknown nodal displacements . Using the Eqn.(3-1> in 
Eqn-(3.7b) one gets 




(3.Sa) 


Here [ B is a matrix of size 3x45 and is given below 


E B ]f = [ B, 


01" ®0S" 




B_ are the sub matrices and is given by, 
Oj 




J 


J 


SH 
a 

^x^ 

c 

dN* 

J 

^X , 


N® 


N. 
) 

R, 


O 


represents column vector containing 
displacements and is given below 


O 


C3.8c) 


unknown nodal 


and 


{ = { *1’ ^E- 

{ 4 - c 


“r “2" “3" ^1' ^2 ) 


}' 


(3.Sd) 


The matrix [ B 3° does not involve the nodal degree of 
freedom. Therefore the first variation of the linear membrane 


strain can be written as 



3C 


4*1 =[«]!'[“•*] 


3. 3. 3. Expressslon tav d 


- V J Mi- 


Rewriting the Eqn. {3.7b) 


{‘F 

JmLi 


1 

2 


du 


3 0 




du. 


dx. 


**1 

r 


dUg 


'*1 


"“3 


dxg 

Jl 

I. 


{3 


^3 

This nonlinear aid plane strain can be written in terns of 
dof as. 


r 1 e e , 

{.} = -[E] [C] {.} 


{3 


Where 


and 


[ G ]« = [G,, G- 


S 3 


E GJ 


aw 


o 


^x 






(3 


(3 


(3.9) 


.10a> 


nodal 


.lOb) 


. lOc ) 


.lOd) 


The equation (3.10b) can be rewritten as 



B o 1 e 

= ; [ m [ •* ] 

Ht. 2 


C3.11a> 


where £ ® ] 
and 


o 

HL. 


C ®01' ®0S 


r B 1 

i- OJ-< 


OJ-lNla 


o 


o 


JT 

^“3 

f 

J 

V 

dx^ 


dN®l 

J 


dx„ 

^ S/ 


^09 3 NL 


o o 


^ ""“a 

iJf 

4* 

^“3] 

J 

d»X . 
1J 



^X_ 

V 




O O 


0 


O 

<3. 11c) 


The first variation of nonlinear mid plane strain can be 
written as 


e 

In [ B matrix, the term is expressed in terms of shape 

functions. After applying variation and rearranging the terms we 
obtain as, 

[ 


Siifastitiitifig tlie relatioii <3-13) in to (3-12) gives 
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r 


r e 


{ « j 

C ® ]1 j 

dd j 

(3.14) 


3. 3. 3 Expression for First, Variation of" Toial Mid Plane Strain s 


The first variation of total mid plane strain can be written 
using relations (3.9) and (3.14) as. 


Where 


.{.}*■ = [Br 

[B]^ = [ [B]«% [B]:/ ] 


(3.15a> 


{3.15b) 


The total mid plane strain is given by. 


{-r = i {^y 


(3. 16) 


3. 3. 4 Expression for d ^ and d-^ e ^ s 


The bending and shear strain at the mid plane of the 
laminate can be written again in terms of nodal displacements as. 


{«} = [B]' {-} 

{‘Y = [■>]'{"}' 


{3.17a) 
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Where 


and 


[ S3 


o 




[ c3" 

== 

[ Cg 

S3 




(3.17b) 



[ D^, Dg 

Q 

1 

1 

1 

0 

0 

0 J 

0 1 




‘'o dx. 


-C 


o dx. 


0 


dx. 


dN 

4 

dx. 




dx 


dN 


dx. 


(3.17c) 




N' 

_4 

R, 


N* 

4 

R, 


O 


^n: 

4 

dx. 


dN 

4 

^X, 


N 


N 


(3.17d) 


The first variation of bending strain and shear strain ii 
given by. 


i-r = 

E 

]' 

1 d * 


E ^ 

]' 

d <5 

Substituting variation 

of 

strain 


(3.18) 


(3.18) in to virtual work equation (3.5) we get. 
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m 

r T 

N -.1 

■T j 

H.l 

yf<5) = 


f e V 

{ d6 } [ Bj" - 

1 

«2 

. +Id6| [C]®. 

1 



- 


V. J 

. «6 


+ 


{ 



dx^ dXg 


C3.193' 


Mhere R contains nodal forces and is 
3.19 can be written after taking out 


derived 



in next section. 


Eg. 


r r 

r N.1 



r r 

1 

T 

1 

VM)= [ B ]® - 

«^2 

- [ c]^. 


J L 



« 6 l 


0 

(3.20> 



Eg. (3.20) represents a set of non linear eguations to be 
solved for a bending problem. 

3.3.5 Expressions For Stress Resultants i 


The expressions for stress resultants ^1 

and can be obtained by substituting the relations (3-16} 
and (3.17a} into Eg- (2.7) we get. 



35 


{ *^ 1 ’ '^6 } 

{ "l- "e- "6}" “ 

{v-^eF = 


[ Oft ] 

[ se] 

I as] 

[ SD ] 

[ QS] 


[c]'{-F 

[B]:{-}*- 


[ «« 3 [ B 31" { -5 }"* 

(3.21a) 

C “ 3 [ Bl] { 5 

(3.21b) 


(3.21c) 


Substituting the expressions 


■for stress resultants ie-. 


relations (3.21a), (3.21b) and (3.21c) in (3.20) we obtain 




J[ 


e’ 1 


[ B ]1 C OA] [ B]; +- [ B]« [ QA] [ B]^' H- 


[B]f [BA] [C]' -H [B]^' [OA] [B]f . 


[ B]1 [ 8B] [ B]“^ . [ B]“^ [ BB] [ C]“ . 

T 1 e 

[ C ] ^ [ SB ] C B ]“' [ C ] ' [ QB ] [ B ]^ ♦ 


[ B] [ «M>] [ C] 


* [ B3'" [ BB] [ D]^] 

{ ■* } '‘* 1 '**E ■ { B “ 


(3.22) 
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3.4, Solution of MonLlnoar Equdlibrlum Eqtmtlons s 


The nonlinear equilibrium equations can be solved by two 
different methods viz. the direct iteration method with 
incremental approach and the incremental newton— raphson method. 

In the direct iteration method load is applied in small 
increments- Some initial solution is guessed (typically, 
6 = equilibrium equations are solved to get new 
solution. This solution is used again to solve the equations. This 
type of iteration is continued till the two solutions converge. 
This method may not converge for highly nonlinear problems. 
Moreover it is computationally expensive. 

In the Newton Raphson method to start with a solution is 
guessed. An improved solution ( d <5^ ) is obtained from 


< 5 .+ d&J 






{3.E3a> 


where K is the tangent stiffness matrix given by. 


E 


' y 

aisy^ 


(3.S3b> 


Eqn- {3.23a> gives nonlinear incremental equilibrium 
equations which gives the linearized approximation to the relation 
between residual forces and the incremental displacements. 

The new improved solution can be found as 
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<5. — S + d <S 

v+4 V i 


(3.S4> 


To improve numerical stability the 
increments. 


load 



is 


applied 


in small 


3.5 Tangent, Stiffness Mat,rlx a 


The tangent stiffness matrix for a element can be written 

as , 


C 





^1 


ss 

M 

^2 

I 

I 

I 

W 

cu 

I 

I 

I 


^E 

M 

m 

m 

m 

Sw 


Sm 

_ « 



SS 

w -1 


Mhere m is the total degree of freedom of the element. 

An element K in Eq.O.SO) can be written as 
Ti.j 
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Tv ,4 


ss 
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^ r 




ss 




SB 


hi 


ss 


ki 




SS 


+ S D 


k- 1 


ki 




ss 


dx^dx„ 
1 2 


C3.26) 


Where 


SR 


SSj 


= o 


It is shown in the Appendix — II that the tangent stiffness 
matrix can be written as 
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e^r Na 1 

e 




1 

[ QS] [ 

1 D ; 

] - 

c <^: 

1 [xl 4 \ 

[ s] 

dx^ dXg 









C3 

.27) 

Using the 

relations 

(3. 

15b) 

in 

Eqn. (3.27) 

the tangent 


stiffness matrix can be written as 


[ x], “ C ><]„ 


Where 


K 


K 


HL 


K 


Linear stiffness matrix 
Initial displacement matrix 
Initial stress matrix 
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[ ><], 
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]“ 1 

‘ QA 3 
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®3f [ 
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E c 
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T 
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3 1 
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[ QB] 

[ B 

]f" [ 
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c]' [ 

GD ] 
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e 

3 1 
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e 

3 

<ix. 
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{3.28b) 
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Jhl< 

f[ b; 
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T 

e 

]: 1 

QA ] 

[ ® iZl- 1 

• -INL 1- 

QA 3 

E ^3:"-^ 

A - 

C c 

T 

e 

3 E 

QB 3 

[ « rZ* 1 
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‘ ® 3*"^ E 

. Jhl. L 

QB 3 

E O" - 
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(3.2Sc) 


K 3 = 

Jo 


T 

m 


A"^ 




«2 1 


(3.28d) 


3. 6 Elemernbal Force Vector i 


The work done by external loads for the bending problea aay 
be independent or combination of these load cases. 

(i) Concentrated loads acting in the direction of the 



corresponding nodal degree of freedom. 

(ii) Uniformly distributed load on the surface. 


4 C 


The total external work done by these forces on an element 
can be written as 


U, 


ER 


W, 


ER 


{ ^ I { ® r 


(3.29) 


R^ is vector containing nodal forces/couples at node j in the 
direction of that nodal degree of freedom. The vector R^ can be 
written as 


and R is 

OJ 


each S 
a 


= {’’or 

given by ^ 

{ s } = { s, . 

is given fay -f O 



(3.30) 


(3.31) 


3.7 (Seometric Stiffness Matrix for the Plate Problem » 


For the linear buckling problem the geometric stiffness 
matrix can be derived by considering the work done by in plane 
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•forces N^and at the mid plane of the laminate. The datum 
configuration is taken as initially compressed configuration. The 
mid plane strains due to bending can be written as 



The work done by these forces is given by 


(3.32) 





"*“3 1 

^ dx.dx 

I ^^<2 J 

<3. 33) 


Eqn. (3.33) can be written as 


W.. 




^“3 I 

T 

1 

H 

dx^ 

1 

- 

2 






- dx„ J 

L 


N. 




N, 


N, 


-1 

^“3 


d»X^ 


1 



•J 

' ^x„ - 


dx. dx_ 
1 2 


(3.34) 


dU dUg 

Expressing the terms and — in terms of shape 

functions and nodal displacements Eqn. (3.34) can be rewritten as 
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1 T T 

- [.5]=[g]^ 

A 

In this expression ,the stress resultants , Ng and are 
constant over the domain and will be equal to applied forces. The 
in plane forces N^and can be expressed as 

= X n^ ; Ng = X ng j ~ ^ (3.36) 

Substituting the relation (3.36) in (3.35) and rewriting the 
expression for external work we obtain 

W® = ^ (3.37) 

Where is the geometric stiffness matrix 

The first variation of the work done by external in plane 
forces can be written by considering the symmetry of as, 

dU* =X [ dfi ]'’ [ K]§ [i ]' (3.38) 

Now the equilibrium equations for a linear buckling problem can be 
written as. 
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or [ [X].' [K ]S ]['*]'= 0 <3.40) 

Eqn. <3.40> represents typical eigen value problem. On 
solvingr the lowest magnitude of eigen value gives the buckling 
load. The eigen vector represents the mode shape. 

3. S Oeome-Lric Sblffness Matrix lor the Shell Problem : 


The geometric stiffness matrix or the shell problem can be 
derived by considering the work due to pre buckling stresses 
Ng, and N^) during the buckling. 

The work done expression can be written as 



2^6 



(3.41) 


Expressing u^ , Ug and Ug in terms of shape functions and the total 
expression will become similarly as Eqn. (3.37) .This Kg matrix is 
substituted in equilibrium equations and solved for the eigen 
value problem. The lowest eigen value is the buckling load. 

The integrals involved in the expressions of ¥'<*5) » 

and K are evaluated using gauss - quardrature. The exact 

& 



44 


integration of the stiffness matriK requires 3x3 gauss point rule. 
But, at the thin plate limit to overcome locking 2x2 gauss 
quadrature is used for shear energy terms, while 3x3 is used for 
non shear terms. In the above section equilibrium equations are 
developed for an element. These equations can be assembled using 
the connectivity at the element interface. Proper boundary 
conditions are applied and the assembled equations are solved. 



CHAPTER 4 


RESULTS AND DISCUSSION 


4. 1 Introduction : 


The equations developed in Chapter 3 are employed for 
numerical studies on some typical plates and shells. The results 
obtained (ie., present linear theory (PLT), present nonlinear 
theory (PNLT) are compared with some experimental, higher order 
theories <HSDT) and classical plate theory (CPT) results. Results 
are presented in the form of graphs and tables. The parameter 
whose effect is examined include number of layers, side to 
thickness ratio (radius to thickness ratio for shell > and 
orthotropy. The effect of different simply-supported boundary 
conditions on bending-stretching coupling is also investigated. 

The finite element used for the discretization of plates 
and shells is a nine noded shell element. In the present work a 
shear correction factor (5/6) has been used. Reduced integration 
has been used for transverse terms in the case of thin plates and 
shells because full integration leads to over stiff results. The 
plate problems are solved by taking the radius of curvature very 
large (typically 10^° ) in strain displacement relations Eqn. 
(2.3a) The boundary conditions used in the present work are shown 
in the Fig. 2.3. 
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For linear problems the matrices K and K , in formulatino 
tangent stiffness matrix are put as zero. For stability problems 
the matrix matrix is formed and Eqn. (3.40) is solved for the 


eigen value. 

For nonlinear bending problems a solution & 
assumed and for the first iteration linear solution 


= { O } is 

is obtained. 


Using the linear solution the matrices K and K are computed and 
solution is found again from the Newton - Raphson iteration. The 
iterations are continued till the two successive solutions j 

3^ -j- 

S converge. The following convergence criterion is used. 



Where i is the iteration counter and k is the nodal 
variable counter. 

Whenever symmetry exists quarter plate model has been 
employed. Unless specified, a 3X3 mesh is used for bending and a 
SX2 is used for buckling problems. The following nondimensional 
deflections and buckling loads are used. 
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Where is calculeled at the center of the plate> N is 
the appiisd critical load in case of linear hucfclin^ prohlesin 

A. 2 Presentation of Results and Discussion t 


The -following types of problems are solved using the 
formulation given in chapter 3. 

4.E.1 Bending Analysis 

(a) Numerical accuracy of the present model with respect to 
the element 

(b) Bending analysis of three layered cross-ply 

(c) Bending analysis of angle-ply 

(d) A comment on the boundary conditions 

(e) Bending of isotropic spherical shell 
<f> Bending of cylindrical shell 

(g) Cross— ply freely supported cylindrical shell under 
sinusoidal load 

4.5.2 Stability Analysis 

(a) Stability of isotropic plate 

(b) Stability of symmetric cross-ply 

(c) Stability of angle-ply 

(d) Stability of isotropic cylindrical shell 

4.2.3 Nonlinear Bending Analysis of Plates 

(a) Nonlinear analysis of orthotropic plate 

(b) Nonlinear analysis of symmetric cross-ply 

(c) Nonlinear analysis of simply supported cross-ply and 


angle— ply 



48 


(d) Nonlinear analysis of clamped corss-ply and angle-ply 

(e) Sandwich plates 

4.2.4 Post Buckling of Isotropic Plate 


4.2.1a Numerical Accuracy of the Present Model with respect to 


the Element : 


To check the numerical accuracy of the present element a 
test problem is solved for which closed form solution (CFS) is 
available C113. 

A simply supported square cross— ply (0/90/90/0) laminate 
with the following properties was considered under sinusoidal 
loading. 


= 25 


V 


12 


0.25 



= 0.5 



= 0.2 , 



(4.1) 


The results of non dimensional ized deflections u^ were 
obtained for the following cases 

(a) Full integration (3X3 gauss rule for bending terms and 
3X3 gauss rule for shear terras) 

(b) Mixed integration (3X3 gauss rule for bending terms and 
2X2 gauss rule for shear terms) 

(c) Reduced integration (2X2 gauss rule for bending and 

shear terms) 

The solutions were obtained using si boundary conditions 
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which are used in Ref C113 for the quarter plate model. 

The results for various thickness ratios were obtained and 
presented in table 4.1. Ue can conclude from the results that 

(1) The nine noded element give virtually same results for full, 
mixed and reduced integration schemes. 

(2) As thickness of the plate is reducing it can be seen that full 
integration gives over stiff results. This is due to the reason 
that the present model being not able to satisfy the Kirchoff’s 
condition as a constraint. It can also be seen that the results 
with mixed integration scheme exactly agree with closed form 
solutions. No significant difference is observed between reduced 
integration and closed form solutions. However the experience of 
the people in this field says that mixed integration is suitable 
for thin plate or thin shell situations. Therefore we used reduced 
integration for the shear terms in these cases. 

(3) The finite element model used here is expected to yield 
correct results. 

4.2.1b Bending analysis of Three Layered Cross-Ply s 


The object of this study is to compare the present theory 
with the closed form solutions of higher order shear deformation 
theory (HSDT) presented by Reddy C163. The higher order theories 
are computationally expensive and involve large number of degrees 
of freedom. 

A simply supported three layered corss-ply (0/90/01 of the 
following material properties subjected to uniformly distributed 



lateral load is taken up for which HSDT results are available 

X 10* psi. Eg = 1.56 X 10* psi, 

^1S ~ ^13 ~ 0.82 X 10 psi, ^23 ~ 0.523 x 10*p5i, 

= 0.24 , a = fa = 10 in. (4.2) 

The results are obtained using s1 faoundry conditions. The 
non dimensional ized deflections Ug are presented for various 
thickness ratios in table 4.2. 

It can be seen from the results presented in table 4.2 that 
the present results are close to the HSDT results for a/h >10 
without much computational effort. The maximum difference is about 
2.7%. The present results underestimate the deflections compared 
to HSDT. It can be seen that the difference reduces as the 
thickness of the plate decreases. The CPT results underestimate 
the deflections and are not reliable at thick plate limit. It can 
be concluded that present theory is adequate in predicting the 
deflections for a/h >10 

4. 2.1c Bending Analysis of Angle-Ply : 


In this section the effect of bending— stretching coupling 
on linear deflections for an angle-ply (-45/45/...) is studied. A 
simply supported angle— ply made of graphite epoxy of the following 
material and geometrical properties subjected to uniformly 
distributed lateral load is taken up. 
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= 40 , Ggg = 0.5 Eg, 

= 0.E5, Eg = 10^ psi 

b =10 in. 


1S 


^13 


(4.3> 


Fig. 4-1 shows the plot of nondimensionalized central 
deflections for various thickness ratios with increasing number of 
layers. The results are obtained using sE boundary conditions in a 
quarter plate. It can be seen that the deflections for two layers 
are higher than the four layer and eight layers. This is due to 
the bending stretching coupling which is maximum for two layers 
and decreases as the number of layers is increased- 

4.2. Id A Comment, on the Etoundary Conditions : 


To study the effect of different in plane boundary 
conditions on bending— stretching problems a two layer angle— ply 
and cross— ply of the material properties given in Eqn. (4.3> 
subjected to uniformly distributed lateral load is taken. 

a b: b = 10 in. h = 0-1 in. C4.4) 

One can use s1 boundary conditions and sE boundary 
conditions as a sya^etry conditions for the quarter plate model 
associated with cross— ply and angle— ply respectively. However one 
cannot use the symmetry conditions si in the quarter model of 
angle-ply and sE in the quarter plate model of cross-ply. The 
effect of different boundary conditions on the linear bending of 

CENT'^''\L. L’?.RARY, 

' ' T., AA'v.ru'R ■' 



40 . 


^23 = 0-5 Eg 



r 


0.25, Eg = 10^ psi 


a = b = 10 in. 


= G^„ 0-6 E„, 
12 13 2 


(4.3) 


Fig. 4.1 shows the plot of nondimensionalized central 
deflections for various thickness ratios with increasing number of 
layers. The results are obtained using s2 boundary conditions in a 
quarter plate. It can be seen that the deflections for two layers 
are higher than the four layer and eight layers. This is due to 
the bending stretching coupling which is maKimum for two layers 
and decreases as the number of layers is increased. 


4. 2. Id A Comaneitt on the Boundary Conditions s 


To study the effect of different in plane boundary 
conditions on bending-stretching problems a two layer angle-ply 
and cross-ply of the material properties given in Eqn. (4.3) 
subjected to uniformly distributed lateral load is taken- 


a = b = 10 in. h = 0.1 in- 


(4.4) 


One can use si boundary conditions and s2 boundary 
conditions as a syaaetry conditions for the quart.r plate aod.l 
associated with cross-ply and angle-ply respectively. However one 
cannot use the syiwsetry conditions si in the quarter aodel of 
angle-ply and s2 in the quarter plate aodel of cross-ply. The 
effect of different boundary conditions on the linear bending of 

L'~RARY 






i 
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cross-ply and angls-ply are showed in the table 4.3. It can be 

seen that for N = 2 layers. thi» •. * 

layers, me difference between the deflections 

obtained from these boundary conditions (S1 and S2> is maKimum. 

But for N = 10 the difference is small. This due to the fact that 

as the number of layers are increased the bending-stretching 

coupling decreases and the effect of different in plane boundary 

conditions is not significant. 

4. 2. le Bending of Isotropic Spherical Shell : 


A simply supported isotropic spherical shell shown in 
Fig. C4. Ha) of the following material and geometrical properties 
is studied under a concentrated load of 100 lbs. 

7 

E = 10 psi, V = 0.3, h = 0.1 in. = 96 in. 

a = b = 32 in. C4-5) 

Due to the existence of the symmetry a 2X2 mesh is used to 
discretise one quadrant of the shell. As the shell is thin (R/h = 
320), mixed integration is used. The linear deflection obtained 
under the concentrated load is 0-03723 in. The mesh is refined to 
3X3 and the central deflection obtained is 0.03886 in- in 
comparison to 0.038661 in. of Rao C373 whose result is based on 
classical shell theory where shear deformation is neglected. 

4-.a.lf Bending of Cylindrical Shell : 

A clamped isotrpic cylindrical shell shown in Fig. (4-Eb) 
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of the following material and geometric properties subjected to 
uniformly distributed lateral load is considered for the study. 

E = 3103 N/mm^, v = 0.3, h = 3.175, a = 254 mm. 

= 1000 X 10~^ N/mm^ e =0.1 rad. (4.6) 

Here a 3X3 mesh is used to descretize one quadrant of the 
shell. Clamped boundary conditions given in Eqn. (2.10) are used. 
The linear central deflection obtained at the center of the shell 
is 0.84297 in. compared to 0.8250 in. of Dutt whose result is 
based on classical shell theory. The Dutt result was given in 
Ref .C393 


4. 2. Ig Cross-ply Freely supported Cylindrical St^ll under 
Sinusoidal Load i 

A cylindrical cross-ply (0/90/90/0) cylindrical shell of 
the material properties given in Eqn (4.1) subjected to 
sinusoidal loading is taken for the study. The geometrical 
properties are 

30 

a/b “ 1-0, a/h = 10 Rg - 10 

R = W SIN SIN -g 

o ,,10 A a 


The following boundary conditions are imposed 
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Ug = Ug = 0 at = 0,a 

Ug = Ui = 0 at = 0,fa (4.S> 

For this special case of boundary conditions , loading closed form 
solutions are available in Ref C383 with the same theory. A 4k4 
element mesh is used to descretize the quarter shell. The 
nondimensional deflections u^ are presented in Table 4.4 for a/h = 
10 and a/h —100- It can be seen that the FEM results are close to 
the closed form solutions. 

4.2.2 Sbablliby Analysis : 

In this section stability of isotropic plates and cross— ply 
plates have been compared with HSDT results. Bending stretching 
coupling on angle-ply laminates has been studied. Results are 
presented for cylindrical shell under axial compressive load. 

4r. 2. 2a Stability of Isotropic Plates : 

The object of this study is to evaluate higher order 
theories with respect to their computational power and accuracy. A 
simply supported isotropic square plate of the following material 
and geometric properties under uniaxial compressive load is 
considered. 

E » lO^psi, V »=s 0.3, a = b = 10 in. (4-9) 

A 2X2 mesh is used to discretize the quadrant of the plate- 
si boundary conditions are used and linear buckling load is 



computed by solving the eigen vaii.*. .. .• 

problem. Nondiroensional 

buckling loads are presented in table 4 s 3 „rf ^ , 

««ie and coiipared with closed 

form solutions of HSDT by Reddy C173 for various thickness ratios. 
It can be seen from the results the present results are very close 
to HSDT results for all thickness ratios. The maximum difference 
is 0.75%. The difference is decreasing as the thickness is 
decreasing. But in the results shown this effect cannot be seen 
over the whole range due to the small error in the eigen value 
calculation of large matrices and locking phenomenon. The CPT 
results over estimate the buckling loads for thick plates (a/h < 
EO). This due to the shear deformation neglected in CPT. 


4 . B. Sb Stability of Symmetric Cross-Plyr : 


A simply supported graphite epoxy symmetric cross-ply 
<0/90/90/0) laminate of the material properties given in Egn. 
(4.3) is considered to study the uniaxiail buckling load. 

The nondimensional uniaxial buckling load was computed 
and presented in table 4,6 with )©DT results C173. si bouniiary 
conditions are used. 

It can also seen from the results that the present results 
are very close to HSDT results. The maximum difference is 1.5%. It 
implies the adequacy of the present t^heory in predicting the 
buckling loads. I«DT in FEH will involve a lot of additional 
computational effort without improvement in the results- From 
Table 4.5 and 4.6 w. can conclude that the difference between the 
present theory and the HH5DT i« slightly increasing 



orthotropy. 


ji.S.2.c: St*billt,y of Anglo-Ply : 


A simply supported square angle-ply t45/-45/.. under 
uniaxial compression load is considered for the study. The 
properties of individual laminate are given in Eqn. (4.3>. 

Solutions are obtained for S2 boundary conditions in a 
quarter plate. Fig. 4,3 displays the variation of nondimensional 
buckling load with fibre orientation and number of layers. It 
is shown II413 that the solution approaches to orthotropic solution 
at N = 6- Hence in the present study N = 6 and N = 2 are 
considered for a/h = 10 and a/h = 100. It can be seen that as the 
number plates increased the buckling load is increasing for the 
same thickness. This is due to the reduction in the bending 
stretching coupling and buckling load is minimum for two layer 
laminates. There is significant difference between the results for 
a/h =* 10 and a/h =» 100. This is due to shear deformation which 
decreases the buckling load. At 45** orientation the buckling loads 
are maximum for all number of layers. 

4.£!.a.<l Buckling of Imotropie Cylindrical Shsll under Axial Load s 

An isotropic homogeneous cylindrical shell with length 11 
in., radius of 36 in., and thickness of 0.125 in. is considered as 
a test problem to check the formulation. Both the edges are simply 
supported. The analytical solution for the same is available with 
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Donnell's theory in Ref C42D. For the finite element modeling it 
is assumed that symmetric conditions prevails at the mid length 
and that the cylinder buckles in 15 ci r cumf erencial full waves as 
shown in the Fig. C4.4a). Consequently a 5.5 in. long and 12*^ wide 
cylindrical segment is considered with symmetrical conditions 
imposed on the three sides. The following boundary conditions are 
used 


^2 

= U3 

= ^2 

= 0 

on 


II 

0 

w 

i:m 



If 

= 

= 0 

on 

^2 

= 0 

1^ 

0 

1 

Ug 

= ^2 

= 0 


on 

^2 

= b 



where coordinate anes pass through left corner of the segment as 
shown in the Fig (4.4b>. 

The eigen value problem of Eqn. (3.42) is solved and 
buckling load is applied as 2721 Ib/in. The result is with a 
difference of 3-5X with the analytical solution (2625 Ib/in.) 
given in Ref.C423. This difference is likely to get reduced as the 
mesh is refined. This is being the reason that as the shell is 
very thin and the shear deformation effect is expected to be 
negligible and the results obtained from the present theory can be 
compared with the Donnell's theory given in Ref. C423 which does 
not take into shear deformation. 

In the present study we could not make use of a routine 
SSPACE given in Ref- C433 which computes specified number of 
lowest eigen values only- Hence we could not go in for a refined 
mesh than 2X2 because of the computational limitation. Laminated 
shells buckles in more than one half wave along the axial 
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was obtained by using both Newton- Raphson and picard iteration. 
Table 4.7 shows the results obtained from Newt on-Raph son and 
picard type of iteration. 

It is evident from Table 4.7 that the number of iterations 
for convergence increases with the load when picard iteration is 
used. However Newton-Raphson converges in three iterations for all 
loads. Hence in the problems investigated below Newton — Raphson 
method has been employed. 

A 2XZ mesh is used in quarter plate. Using Newton-Raphson 
results were obtained with S1 and S3 boundary conditions. Fig. 
(4.5) displays the results obtained from present linear results 
(PLT), present nonlinear results (PNLT), classical linear results 
(CLT) and classical nonlinear results (Cl'fl-T) given in Ref. C183. 
It can be seen that the PLT and CLT results are close. This is due 
to the fact that as the plate is very thin (a/h = 87> the effect 
of shear deformation is expected to be negligible- It can be seen 
that compared to CNLT results the PNLT results obtained with S1 
boundary conditions are in close agreement with the experimental 
results. This is due to the reason that the boundary conditions 
used in Ref. C183 to obtain the ChH-T solution are = Ug = 
u„ = O and M„ » O at the simply supported edge which may not 
correspond to the the experimental conditions. The Pr«-T results 
obtained with S3 boundary conditions are away from the 
experimental results. The difference between S1 and S3 boundary 
conditions are due to the bending stretching introduced through 
the nonlinear term. 

It can be concluded that the practical simply-supported 



boundary conditions for the orthotropic plate can be in between SI 
and S3 boundary conditions. 

4.2.3b Nonlinear analysis of Cross-Ply s 


A clamped symmetric cross-ply (0/90/90/0) of the following 
material and geometric properties used to study the nonlinear 
analysis under uniformly distributed load. For this data CLT, CNLT 
and experimental results are avalible in CISIl. 


b = 12 in, 


E, = 1,828 X 10 
1 


h = 0.096 in. 

Eg = 1.8315 X 10^ 


V 


12 = “13 =“23 = 


<4.12) 


The CRT, CNLT, PLT, PNLT and experimental results are 
presented in Fig. (4.4) with load intensity. In this case there is 
a difference between CUT and PLT reults even though the plate is 

.. ies tho Ref C183 the following boundary 

thin. This is because in the Ref. 

^ sst th- clamoed edge in the linear and 
conditions are used at the damp u 

nonlinear analysis 


du. 


(4.13) 


= “2 * “3 


Sx. 


which over estimates the deflections. But in the present FEM model 
it is possible to specify the correct boundary conditions which 

■ Pnn (2 10) The PNLT results are lower than the 

are given .m Egn • 

*4. .rsm-e-Oe-t t T SO 4 . Th 6 T 005 OO Of 
experimental results which shows e 



61 


difference in the results is that the clamped boundary conditions 
are not satisfactorily imposed during the experimentation as 
stated in the Ref.C18j|. The CNLT results shows altogether 
different trend because of the reason mentioned above regarding 
boundary conditions. These two examples led to the confidence in 
the present model in predicting the accurate nonlinear response of 
the composite plates. 

< 4 . 2. 3c Nonlinear- Analysis of Singly Supported Cross— Ply and 
Angle-Ply i 


A simply supported square plate of the following material 
properties is taken for the study. 



, Eg = 7.031 x 10^ N/cffl^ , 


23 


= .2 


*^12 


243.8 cm. 


h = .635 



.5 


.25 


(4.15) 


As the bending— stretching coupling is maximum for two layer 
laminates a two layer angle-ply and cross-ply subjected to 
uniformly distributed lateral load is considered. The load 
deflections curves are shown in Fig- (4.7) and Fig. (4.8). For the 
cross-ply SI and S3 boundary conditions are used. It can be seen 
that S3 boundary conditions are stiffening type than SI boundary 
conditions. The same trend is shown for the orthotropic plate in 
Fig. (4.4). For the angle-ply S2 and S3 boundary conditions are 
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used- Here also S3 boundary conditions show the same trend with S2 
boundary condition. 


4,. 4.. 3d Nonlinear Analysis of Clanged Cross-ply and Angle-Ply s 


A clamped cross— ply and angle— ply of the material 
properties shown in Eqn.(4.3> is used to study the nonlinear 
behavior . 

Fig. (4.9) and Fig. (4. 10) shows the load deflection curves 
for two layer cross-ply and angle-ply. It can be seen that the two 
layer cross-ply exhibits more nonlinearity than angle-ply. It can 
also seen that for cross-ply as the number of layers are increased 
from two layers to eight layers the nonlinearity is decreasing. 
This is due to the reduction in bending-stretching coupling with 
the number of layers. 

A.3.3& Sandwich Plates i 


The object of this study is to validate the present theory 
to sandwich plates by comparing with the available experimental 
results. A simply supported square sandwich laminate with 
aluminium faces and balsa core subjected to uniform transverse 
pressure is studied with the following face properties. 

« 10 X 10^, = 0-32 in. 

A 2x2 mesh is used to idealize one quarter plate and S3 
boundary condition are used. The results are obtained from the 

present theory and experimental results are presented in Table 4.8 



and table 4.9 for varying spans, core thickness, core shear 
rigidity and load intensities. The present nonlinear results 
compare satisfactorily with the experimental results. When the 
shear rigidity of the core high it can be concluded that 
neglecting the face shear rigidity yields better results than 
considering the shear rigidity of the face plate. 


4. S. 4. Post buckling of* Isotropic Plate : 


A simply supported isotropic square plate of the following 
material properties is chosen for the study 


E = 30 X lO^Psi, a = b = 10 in. h = 0.1 in. 

In this problem following type of boundary conditions were 


used 




nx. 





SIN 

I 

cos 

“SE" 




ftH. 


nx„ 



COS 

1 

si 

SIN 

b ^2 

(4.165 



fXK. 


nx„ 


= 

COS 

m 

COS 

-HT 



Here the origin is chosen at the centre of the plate- 2a and 2b 
are the dimensions of the plate respectively. 

The boundary conditions can be rewritten as 


**2 * ^2 
Ug = + « b 


on ail edges, 
on ^ 
on x„ = 


(4.175 


+ b 



After the application 


of these boundary conditions the 
problem becomes the minimization of strain energy as the variation 
of potential energy is zero. 

From Ref. C42ZI the critical value of ■£ where a buckled 
configuration exists for h> = 0.3 is 


« ^ = 0 , 632 — = 


( 4 . 18 ) 


ie., if « < «^^the plate will be in compressed configuration only 


Let n = — - 

cr 

A small lateral load 21bs (0.18X of the critical load) is 
applied at the centre of the plate to initiate buckling. The 
buckled configuration serves as starting vector for the iteration. 

The central deflection computed for various values of n 
iteratively are presented with analytical Gelerkian solution 
obtained in Ref. C423 in Fig- 4.11. It can be seen that the 
present FEM results and Gelerkian solution are agreeing. 















Nondimensij^nalized deflections for symmetric cross-ply ( 0 / 90 / 90 / 0 ) 
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r^Lle 4-‘7 Results obtained for Newton-Raphson and 
and Picard Iteration for the Orthotropic 
Plate 


Load 

No . of 

iterations 

Picard 

Newton-Raphson 


1 

4.5365 E-02 

4.5365 E-02 


2 

4.3818 E-02 

4.4142 E-02 

0.2 





3 

4.3928 E-02 

4.3917 E-02 


4 


4.3928 E-02 


1 

8.7836 E-02 

8.5241 E-02 


2 

8.0171 E-02 

8.2001 E-02 

0.4 

3 

8.1796 E-02 

8.1515 E-02 


4 

8.1467 E-02 

8.1522 E-02 


5 

8.1578 E-02 

8.1522 E-02 


1 

0.1223 

0.1154 


2 

0.1085 

0.1125 

o.fc 

3 

0.1135 

0.1122 


4 

0.1117 



5 

0.1123 



6 

0.1121 



1 

0.1497 

0.1399 


2 

0.1316 

0.1377 


3 

0.1405 

0.1376 


4 

0.1362 


0.8 

5 

0.1383 



6 

0.1372 



7 

0.1377 



8 

0.1375 



1 

0.1720 

0.1608 


2 

0.1514 

0.1592 


3 

0.1638 

i 0.1591 


4 

0.1563 



5 

0.1607 


1.0 

6 

0.1581 



7 

0.1597 



8 

0.1593 




0.1589 














Table 4.8 Sandwich plates comparison with experimental results ( G^=12670ps i ) 
A - G.not considered B - G. considered 
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Table 4.9 Sandwich plates comparison with experimental results (G =5840psi) 
A - G.nat considered B - G.considered ^ 
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CHAPTER 5 


CONCLUSIONS 


5.1 Conclusions : 

A first order shear deformation theory for the nonlinear 
bending, post buckling of plates, linear bending and buckling of 
shells is presented- On the basis of the limited study carried 
pout the following conclusions can be drawn. 

1) It is clear from the Table 4.1 that the type of 
integration scheme has a little effect on the nine noded 
lagrangian element. It is also evident that one has to use reduced 
integration for the shear terms in the thin plate region to 
circumvent locking. 

2) It can be seen from the Table 4.2 that FSDT results are 
not reliable for a/h < 10 for the bending problems. For a/h > 10 
there is no significant difference between FSDT and HSDT results. 

3) The FSDT results for the uniaxial buckling load are 
studied for three layered symmetric cross-ply and is given in 
Table 4.6. For this type of laminate studied it can be seen that 
present results are very close to HSDT results for all thickness 
ratios. It can be concluded that the present theory is adequate to 
determine the lowest buckling load. The classical plate theory is 
not reliable for a/h < 20. 

4> It can be seen from the Table 4.7 that NewtonH?aphson 
method is more efficient than the picard type of iteration through 
the former involves computation of derivative of the stiffness 
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matrix. 

5> As seen from the Fig. 4.4 it can he concluded that the 
present nonlinear results are close to the experimental results 
with S1 boundary condition than S3 boundary condition- The 
difference in results for S1 and S3 boundary conditions is due to 
the bending-stretching coupling introduced through the 
nonlinearity. The linear results over estimate the deflections and 
are not reliable where the deflections are of the order of the 
thickness. It can noted that the experimental boundary conditions 
may lie in between S1 and S3 boundary conditions for an 
orthotropic plate. 

6) Due the bending stretching coupling the different in 
plane boundary conditions (S1,SS,S3) give different results. Of 
all the boundary conditions S3 is stiffening type. This can be 
seen from Fig. 4.7 and 4.8. 

7> From Fig. 4.9 and 4.10 it can seen that cross— ply 
exhibits higher nonlinearity than angle-ply laminates. It can also 
noted that for the cross— ply as the number of layers is increased 
the nonlinearity is decreasing. 

5. 2 Reconumendatlons for the Future Work : 

Some of the future work that can be carried out as an 
extension or improvement over the present investigation are 
suggested below. 

The buckling formulation and program that is developed can 
be utilized to analyze the stability of sandwich plates. The 
program developed here also caters for biaxial and shear buckling 
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of laminated composite plates. This could be investigated for 
future work. 

With the little effort this work can be extended to 
nonlinear analysis by taking all the nonlinear terms in the strain 
displacement relation as 


du^ 

4 - 

1 
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I 

r*f 

.U, , 

3u^ .2 
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1 

J t 

1 

dx^ J _ 


e 

There will be additional matrices such as give in 
Eqn.(3. ). This becomes computationally expensive. One can 
investigate whether that effort is required in engineering design. 

The post buckling problem done in the present work is 
limited to some special boundary conditions. By adding an 
additional external potential term one can solve it for all 
boundary conditions. 
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APPENDIX I 


1 Int-roduct-ion s 


A composite laminate is composed by bonding together a 
number of lamina. In each lamina material aKes may not coincide 
with the reference axes. The material constants of a lamina which 
are functions of engineering constants and orientation of fibre « 
are given in the following section. 

2 Description : 

The stress strain relations for a lamina, with material 
axes coinciding with reference axes and assuming the transverse 
stress is negligible, are given by 
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For an arbitrarily oriented laminated lamina after carrying 
out necessary transformation stress-strain relations are given by. 
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In which C= COS(«) and S = SINC«) and « is the angle between 

reference axes and principal, orthotropic axes 



APPE^©IX II 


1. Derivation of Tarigent Stiffness Matrix : 


Substituting the Eqn. C3.20) , in to Eqn. (3.S5) enables 
the tangent stiffness matriK to be written as 
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Where n denotes total degree of freedom. 

Since the stress resultants are linear function of 
the variation of stress resultants can be writen 
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Eqn.(2.7) as 
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expression for the tangent stiffness matrix as 



[ D]® [ QS] [ D]® dA 

(A E.8) 

Considering the third integral in Eqn. (A S.8) and noting 
that [ B is not a function of <5 gives, upon mult i pi cat ion of 

li(^ N ) and E 
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The differential with respect to S in Eqn. CA 2.13) are the 
coefficients of [ G ] matrix of Eqn- (3.10d) typically for a 
degree of freedom i. 
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Now J'M(N) can he written as 
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Consequently Eqn. (A 2-17) may be reaaranged to give the 


expressions for K shown in Eqn C3.27a) 
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